 /* partInfo.c

   Usage: ./partInfo [-h help] <.graph file> <.part file>
   Give info about a partition. Take
   - a graph file and
   - a partitioning file (metis or pagrid style)
   - optionally a vmap file (pagrid input's processor grid)

   Basile Clout, November 2007
   */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <ctype.h>
#include <unistd.h>

#define LINE 512
#define MAX_DOMAINS 128
#define MAX_WEIGHT 999999999
#define G 10
#define LW_CEIL 1000		/* Max latency. Taken from pagridL (grid_func.c) */

/*#define __DEBUG__*/


int gwBdw(int **sn, int **wbdw, int i, int j);
int gwLat(int **sn, int **wbdw, int i, int j);


int usage(){
  printf("Usage:\n");
  printf(" partInfo [-h help] <.graph file> <.part file>\n");
  printf("  - graph file must be an adjacency graph (such as created with mesh2graph.py)\n");
  printf("  - part file is a metis-style partition file\n\n");
  printf(" partInfo [-h help] [-p] [-g G pagrid parameter] <.graph file> <.part file> <.vmap file>\n");
  printf("  - graph file is an adjacency graph (such as created with mesh2graph.py)\n");
  printf("  - part file is a metis-style partition file\n");
  printf("  - vmap file is a vmap pagrid input file\n\n");
  
  return 1;
}
  

int main(int argc, char *argv[]){


  FILE *graph;
  FILE *part;
  FILE *vmap;
  int *ap;			/* partition's array */
  int *adj;			/* Complete array */
  int **pta;			/* Pointer on nodes in complete array ap */
  int *nn;			/* number of edges for each node n */
  int n, e, p;			/* Number of nodes, edges, and partitions */
  char line[LINE];
  int g=G;
  char *myline = malloc(LINE*sizeof(char));
  int *domains = malloc(MAX_DOMAINS*sizeof(int)); /* Contains the list of domains */
  int *cep = malloc(MAX_DOMAINS*sizeof(int)); /* Contains the cut edges from domain i */
  int *np = malloc(MAX_DOMAINS*sizeof(int));  /* Number of nodes in domain i */
  double *eet = malloc(MAX_DOMAINS*sizeof(double)); /* Estimated Execution Time of each processor */
  double *stats = malloc(MAX_DOMAINS*3*sizeof(double));				    /* EET stats for the partition: stats[0]=proc, stats[1]=bdw, stats[2]=lat */
  int **lnp;				      /* List of nodes per processor */
  char *buf[128];
  int m;
  int i, j, c;
  int s;
  double s1, s2;
  int r, v;
  double t;
  int ce=0;
  double maxeet=0;
  int maxk=0;
  double rref;  
  int h = 0;
  int k;			/* One more iterator */
  int vp, ve;			/* Vmap number of processors and edges */
  int p1, p2;			/* Vmap ref node 1, ref node 2*/
  int *vpp, *vadj;		/* Number of directly connected processors list, adjacent processors */
  int *vpbdw, *vplat; /* Arrays containing cpu, bandwidth and latency weights */
  double *vpcpu;
  int maxlat=0;			/* Max cpu. Use to scale back CPU values */
  double rlat;
  double maxcpu=0;
  int **vpta, **vptb, **vptl;			/* Node Pointer in adjacency array, bandwidthand latency arrays (use with nn) */

  int **spaths, **sn;			/* p*p shortest path matrix, and next hop */
  int **wlat, **wbdw;		/* Bandwidth and latency weights between each pair of nodes */
  double **rpr;			/* R matrix (Bandwidth ratios) */
  double **spr;			/* S matrix (Latency ratios) */

  int **cevp;
  

  /* Getopt variables */
  int  gc;
  opterr=0;

  while((gc=getopt(argc, argv, "hpg:"))!=-1)
    switch(gc){
    case 'h':
      usage();
      exit(0);
      break;
    case 'p':
      h=1;
      break;
    case 'g':
      g=(int) atoi(optarg);
      break;
    case '?':
      if(optopt=='g')
	fprintf(stderr, "Option -%c requires an argument.\n", optopt);
      else if (isprint (optopt))
	fprintf(stderr, "Unknown option `-%c'.\n", optopt);
      else
	fprintf(stderr, "Unknown option character `\\x%x'.\n", optopt);
      return 1;
    }

    if (h){
      if (argc-optind != 3){
	printf("[heterogeneous grid] Please enter a .graph, .part and .vmap file\n");
	usage();
	exit(1);
      }
    }else{
      if(argc-optind != 2){
	printf("[homogeneous grid] Please enter a .graph and .part file\n");
	usage();
	exit(1);
      }
    }
  
  if((graph = fopen(argv[optind], "r")) == NULL){
    perror("Impossible to open .graph file");
    exit(1);
  }
  if((part = fopen(argv[optind+1], "r")) == NULL){
    perror("Impossible to open .part file");
    exit(1);
  }
  if(h){
    if((vmap = fopen(argv[optind+2], "r")) == NULL){
      perror("Impossible to open .vmap file");
      exit(1);
    }
  }

  /* Parse graph file */
  if(fgets(line, LINE, graph) != NULL){
    sscanf(line, "%d %d", &n, &e);
  }

  /* malloc require place for adjacency array */
  adj = (int *) malloc(2*e*sizeof(int)); /* Each edge appears two times in the array ... use more memory, but is faster */
  nn = (int *) malloc(n*sizeof(int));
  pta = (int **) malloc(n*sizeof(int *));
  pta[0]=adj;

  /* Fill the adjacency array */
  i=j=0;
  s=0;
  while ((fgets(line, LINE, graph))!=NULL){
    myline = line;
    c=0;
    while ((*buf = strsep(&myline, " \t\n")) != NULL){
      if (**buf != '\0'){
#ifdef __DEBUG__ 
	printf("'%d' ", atoi(*buf));
#endif
	adj[j++]=(int)atoi(*buf)-1; 
	c++;
      }
    }
    nn[i]=c;
    s+=c;
    i++;
    if(i<n)
      pta[i] = adj+s;
#ifdef __DEBUG__
    printf("\n");
#endif
  }

#ifdef __DEBUG__
  /* Print the adjacent nodes of each node */
  c=0;
  for(i=0;i<n;i++){
    printf("\nNode %d:\n", i+1);
    for(j=0;j<nn[i];j++){
      printf("%d(%d) ", adj[c+j]+1, pta[i][j]+1);/* arr start at 0, metis nodes at 1... */
    }
    c+=nn[i];
  }
  printf("\n");
#endif
  
  /* Parse partition file */
  if((ap = (int *) malloc(n*sizeof(int)))==NULL){
      perror("Impossibe to malloc processor space");
    }

  p=0;
  for(i=0;(fscanf(part, "%d", &ap[i])) != EOF; i++){
    m=0;
    for(j=0;j<p;j++){		/* Check if the current domain is the new domain! */
      if(ap[i]==domains[j]){
	np[ap[i]]++;
	m=1;
	break;
      }
    }
    if(m==0){			/* Update the array of domains and the number of different domains p*/
      domains[p++]=ap[i];
      np[ap[i]]=1;
    }
#ifdef __DEBUG__
    printf("domain = %d/%d\n", ap[i], p);
#endif
  }

  if(p>MAX_DOMAINS){
    printf("Too many domains (%d domains, max is %d)\n", p, MAX_DOMAINS);
    exit(1);
  }
  
  /* Build an array for each node with the list of its vertices */
  lnp = (int **) malloc(p*sizeof(int*));
  for(i=0;i<p;i++){
    lnp[i] = (int *) malloc(np[i]*sizeof(int));
  }
  for(i=0;i<p;i++){
    c=0;
    for(j=0;j<n;j++){
      if(ap[j]==i){
	lnp[i][c++]=j;
      }
    }
  }

#ifdef __DEBUG__
  for(i=0;i<p;i++){
    printf("In domain %d:\n", i);
    for(j=0;j<np[i];j++){
      printf("%d ", lnp[i][j]+1);
    }
    printf("\n");
  }
#endif
 

  
#ifdef __DEBUG__
  printf("Cut edges:\n");
#endif
  ce=m=0;
  /* Determine the number of cut edges c*/
  
  /* cut edges for each partition (adjacency) */
  cevp = (int **) malloc(p*sizeof(int*));
  for(i=0;i<p;i++)
    cevp[i]=(int *) malloc(p*sizeof(int));

  /* Initialization */
  for(i=0;i<p;i++){
    for(j=0;j<p;j++)
      cevp[i][j]=0;
  }
  
  for(i=0;i<n;i++){
    for(j=0;j<nn[i];j++){
      if(ap[i]!=ap[adj[m+j]]){
#ifdef __DEBUG__
	printf("(%d, %d)\n", i+1, adj[m+j]+1); /* Print in the metis way: start from 1 */
#endif
	cep[ap[i]]++, cep[ap[adj[m+j]]]++;
	ce++;
	cevp[ap[i]][ap[adj[m+j]]]++;
      }
    }
    m+=nn[i];
  }
  
  
  /* Determine the EET for each processor */
  if(h){

    /* Parse vmap file */
    i=j=0;
    s=0;

    while((fgets(line, LINE, vmap)) != NULL){

      if(line[0]!='%'){
	if(!j){
	  sscanf(line, "%d %d", &vp, &ve);
	  j++;
	}else if (j==1){
	  sscanf(line, "%d %d %lf", &p1, &p2, &rref);
	  j++;
#ifdef __DEBUG__
	  printf("vp=%d, ve=%d\n", vp, ve);
	  printf("n1=%d, n2=%d, rref = %.5f\n", p1, p2, rref);
#endif
	  vpta = (int **) malloc(vp*sizeof(int *));
	  vptb = (int **) malloc(vp*sizeof(int *));
	  vptl = (int **) malloc(vp*sizeof(int *));
	  vpp = (int *) malloc(vp*sizeof(int));
	  vadj = (int *) malloc(2*ve*sizeof(int));
	  vpcpu = (double *) malloc(vp*sizeof(double));
	  vpbdw = (int *) malloc(2*ve*sizeof(int));
	  vplat = (int *) malloc(2*ve*sizeof(int));
	  
	  if(vpta==NULL || vpp==NULL || vadj==NULL ||
	     vpcpu==NULL || vpbdw==NULL || vplat==NULL
	     || vptb==NULL || vptl==NULL ){
	    printf("Impossible to malloc space.\n");
	    exit(1);
	  }
	  vpta[0]=vadj;
	  vptb[0]=vpbdw;
	  vptl[0]=vplat;
	  
	}else{

	  myline = line;
	  c=k=0;
	  
	  while((*buf = strsep(&myline, " \t\n")) != NULL){
	    if(**buf != '\0'){
	      switch(k){
	      case 0:
		vpcpu[i] = atof(*buf);
		if (vpcpu[i]>maxcpu)
		  maxcpu=vpcpu[i];
		k++;
		break;
	      case 1:
		vadj[s+c] = atoi(*buf);
		k++;
		break;
	      case 2:
		vpbdw[s+c] = atoi(*buf);
		k++;
		break;
	      case 3:
		vplat[s+c] = atoi(*buf)/g;
		k=1;
		c++;
		break;
	      }
	    }
	  }
	  vpp[i]=c;
	  s+=c;
	  i++;
	  if(i<p){
	    vpta[i] = vadj+s;
	    vptb[i] = vpbdw+s;
	    vptl[i] = vplat+s;
	  }
	}
      }
    }
  

    if(p!=vp){
      printf("Error: processor graph does not correspond to partition file (p=%d, vp=%d).\n",p, vp);
      exit(1);
    }
    
  /* Scale back all cpu weights (max=1, instead of the conf file hypothesis min=1) */
    for(i=0;i<p;i++)
      vpcpu[i]=maxcpu/vpcpu[i];
    
#ifdef __DEBUG__
    printf("\n");
    s=0;
    for(i=0;i<vp;i++){
      printf("domain %d: cpu=%f vpp=%d\n", i, vpcpu[i], vpp[i]);
      for(j=0;j<vpp[i]; j++){
	printf("-->%d: %d %d\t\t", vadj[s+j], vpbdw[s+j], vplat[s+j]);
	printf("(with pointers:-->%d: %d %d)\n", vpta[i][j], vptb[i][j], vptl[i][j]);
      }
      s+=vpp[i];
    }
#endif
    
    /* Shortest path matrix */

    /* Initialization */
    spaths = (int **) malloc(p*sizeof(int *));
    sn = (int **) malloc(p*sizeof(int *));
    wbdw = (int **) malloc(p*sizeof(int *));
    wlat = (int **) malloc(p*sizeof(int *));
    for(i=0;i<p;i++){
      spaths[i] = (int *) malloc(p*sizeof(int));
      sn[i] = (int *) malloc(p*sizeof(int)); 
      wbdw[i] = (int *) malloc(p*sizeof(int));
      wlat[i] = (int *) malloc(p*sizeof(int));
    }
    
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	if(i==j){
	  spaths[i][j]=0;
	  wbdw[i][j]=0;
	  wlat[i][j]=0;
	}else{
	  spaths[i][j]=MAX_WEIGHT;
	  wbdw[i][j]=MAX_WEIGHT;
	  wlat[i][j]=MAX_WEIGHT;
	}
	sn[i][j]=-1;
      }
    }
    
    for(i=0;i<p;i++){
      for(j=0;j<vpp[i];j++){
	wlat[i][vpta[i][j]]=wlat[vpta[i][j]][i]=vptl[i][j];
	wbdw[i][vpta[i][j]]=wbdw[vpta[i][j]][i]=vptb[i][j];
	spaths[i][vpta[i][j]]=1;
      }
    }

    /* Shortest paths */
    for(k=0;k<p;k++){
      for(i=0;i<p;i++){		/* Can skip first step of Dantzig: adjacency is always direct (all weights are 1) */
	for(j=0;j<p;j++){
	  if((c=spaths[i][k]+spaths[k][j])<spaths[i][j]){
	  spaths[i][j]=c;
	  sn[i][j]=k;
	  }
	}
      }
    }
    
    /* bandwidth and latency weights */
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	if (sn[i][j]==-1)		/* Directly connected (or no paths ...): we already know the weight */
	  continue;
	if(wbdw[i][j]==MAX_WEIGHT) /* Lets find the corresponding weight */
	  wbdw[i][j]=wbdw[j][i]=gwBdw(sn, wbdw, i, j);
	if(wlat[i][j]==MAX_WEIGHT)
	  wlat[i][j]=wlat[j][i]=gwLat(sn, wlat, i, j);
      }
    }

    
    /* Generate matrices R and S */
    rpr = (double **) malloc(p*sizeof(double*));
    spr = (double **) malloc(p*sizeof(double*));
    for(i=0;i<p;i++){
      rpr[i]=(double *)malloc(p*sizeof(double));
      spr[i]=(double *)malloc(p*sizeof(double));
    }
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	rpr[i][j]=(double)rref*(vpcpu[p1]/wbdw[p1][p2])*((double)wbdw[i][j]/vpcpu[i]);
	spr[i][j]=(double)rref*(vpcpu[p1]/wbdw[p1][p2])*((double)wlat[i][j]/vpcpu[i]);
	if(spr[i][j]>maxlat)
	  maxlat=spr[i][j];
      }
    }

    
    if(maxlat>LW_CEIL)
      rlat = (double)maxlat/LW_CEIL;
    else
      rlat=1;
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	spr[i][j]=(double)(spr[i][j]/rlat);
      }
    }
    
    
#ifdef __DEBUG__
    /* Debug stuff */
    printf("Shortest paths\n");
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	printf("%d ", spaths[i][j]);
      }
      printf("\n");
    }
    printf("\n");
    
    printf("Next hop\n");
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	printf("%d ", sn[i][j]);
      }
      printf("\n");
    }
    printf("\n");
  
    printf("Latency:\n");
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	printf("%d ", (int)spr[i][j]);
      }
      printf("\n");
    }
    printf("\n");

    printf("Bandwidth\n");
    for(i=0;i<p;i++){
      for(j=0;j<p;j++){
	printf("%f ", rpr[i][j]);
      }
      printf("\n");
    }
    printf("\n");

    printf("CPU\n");
    for(i=0;i<p;i++){
      printf("%d: %f\n", i, vpcpu[i]);
    }
    printf("\n");    
#endif

    for(k=0;k<p;k++){
      t=np[k];			/* Weight each node=1 */
      stats[3*k+0]=np[k]*vpcpu[k];
      s1=0.0;
      for(v=0;v<np[k];v++){
	for(r=0;r<p;r++){
	  j=0;
	  if(r!=k){
	    for(c=0;c<nn[lnp[k][v]];c++){
	      if(ap[pta[lnp[k][v]][c]]==r){
		j++;		/* We suppose here that weight of all edges is 1 ... */
	      }
	    }
	    s1+=j*rpr[k][r];
	  }
	}
      }
      t+=s1;
      stats[3*k+1]=s1*vpcpu[k];
      c=0;
      s2=0.0;
      for(r=0;r<p;r++){
	if(r!=k){
	  c=0;
	  for(v=0;v<np[k];v++){
	    for(j=0;j<nn[lnp[k][v]];j++){
	      if(ap[pta[lnp[k][v]][j]]==r){
		c=1;
		break;
	      }
	    }
	  }
	  if(c==1){
	    s2+=(int)spr[k][r];
#ifdef __DEBUG__
	    printf("%d: exec=%f. Add latency %d to (%d, %d)\n", k, t, (int)spr[k][r], k, r);
#endif
	  }
	}	
      }
      t+=s2;
      stats[3*k+2]=s2*vpcpu[k];
      t*=vpcpu[k];
      if(t>maxeet){
	maxeet=t;
	maxk=k;
      }
      eet[k]=t;
    }
  }

  printf("Partition statistics. \n- graph \'%s\' \n- partition \'%s\'\n", argv[optind], argv[optind+1]);
  if(h)
    printf("- vmap '%s'\n", argv[optind+2]);
  printf("# vertices:\t%d\n", n);
  printf("# edges:\t%d\n", e);
  printf("# domains:\t%d\n", p);
  if(h){
    printf("# domain edges:\t%d\n", ve);
    printf("G:\t\t%d\n", g);
    printf("Max EET:\t%.2f", maxeet);
    printf("\t(CPU: %.2f, Bandwidth: %.2f, Latency: %.2f, Ratio link/cpu: %.2f)\n", stats[3*maxk+0], stats[3*maxk+1], stats[3*maxk+2], (stats[3*maxk+1]+stats[3*maxk+2])/stats[3*maxk]);
  }
  printf("# CE:\t\t%d\n", ce/2);
  printf("Per domain statistics:\n");
  if(h){
    printf("\tDomain\t#nodes\tEET\t#CE\tAdjacents\n");
    for(i=0;i<p;i++){
      printf("\t%d\t%d\t%.2f %.2f %.2f %.2f %.2f\t%d\t", i, np[i], eet[i], stats[3*i], stats[3*i+1], stats[3*i+2], (stats[3*i+1]+stats[3*i+2])/stats[3*i], cep[i]/2);
      for(j=0;j<p;j++){
	if(cevp[i][j]>0)
	  printf("%d(%d), ", j, cevp[i][j]);
      }
      printf("\n");
    }
  }else{
    printf("\tDomain\t#nodes\t#CE\tAdjacents\n");
    for(i=0;i<p;i++){
      printf("\t%d\t%d\t%d\t", i, np[i], cep[i]/2);
      for(j=0;j<p;j++){
	if(cevp[i][j]>0)
	  printf("%d(%d), ", j, cevp[i][j]);
      }
      printf("\n");
    }
  }
    
  
  exit(0);
}
  

/* Get bandwidth weight */
int gwBdw(int **sn, int **wbdw, int i, int j){
  int bik, bkj, k;
  bik=bkj=0;
  if(sn[i][j]==-1){		/* Lowest level: direct adjacency */
    return wbdw[i][j];
  }else{
    k=sn[i][j];
    bik=gwBdw(sn, wbdw, i, k);
    bkj=gwBdw(sn, wbdw, k, j);
    if(bik>=bkj)
      return bik;
    else
      return bkj;
  }
}

/* Get latency weight */
int gwLat(int **sn, int **wlat, int i, int j){
  int lik, lkj, k;
  if(sn[i][j]==-1){		/* Lowest level: direct adjacency */
    return wlat[i][j];
  }else{
    lik=lkj=0;
    k=sn[i][j];
    lik=gwLat(sn, wlat, i, k);
    lkj=gwBdw(sn, wlat, k, j);
    return lik+lkj;
  }
}

  
